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1.  INTRODUCTION  AND  OVERVIEW 


The  clarification  of  the  complex  physical  phenomena  occurring  in  a 
snowstorm  is  of  great  interest,  both  scientifically  and  from  a  practical  point 
of  view.  On  the  one  hand,  most  of  the  snow  cloud  microphysics  is  repeated  in 
the  upper  levels  of  a  raincloud,  and  is  therefore  of  very  general 
meteorological  interest.  On  the  other  hand,  the  strong  impact  of  falling  snow 
on  visibility  and  transmission  of  various  types  of  radiation,  and  the 
variability  of  these  effects  depending  on  intensity  and  type  of  snow1  indicate 
the  need  for  some  understanding  of  the  underlying  mechanisms.  In  order  to 
develop  a  falling  snow  model  describing  the  adverse  weather  effects  of  snow  on 
transmission,  it  is  necessary  to  incorporate  optical,  physical,  and 
meterological  characteristics.  We  report  on  development  of  a  model  which 
incorporates  these  aspects  of  snow  characterization. 

In  Sections  2  to  5,  we  report  on  the  development  of  a  numerical  model  of 
snowfall  aimed  at  obtaining  such  an  understanding.  Specifically,  the  model 
accepts  as  inputs  such  measurable  quantities  as  temperature  and  dew  point 
vertical  profiles,  as  well  as  other  quantities,  often  not  measured,  such  as 
updraft  velocity,  cloud  droplet  size  and  ice  nucleus  concentration.  Based  on 
these  data,  the  model  then  predicts  the  time  evolution  of  the  cloud,  and, 
eventually,  the  rate  of  precipitation  and  the  type  of  snow  particles  arriving 
at  the  ground.  This  then  provides  some  degree  of  closure,  by  allowing 
comparisons  to  such  readily  observable  quantities  as  equivalent  water 
precipitation  rate  and  particle  flux  at  the  ground. 

The  model  can  display  the  sensitivity  of  the  snowstorm  to  specific 
changes  in  the  data,  and  can,  therefore,  serve  as  a  guide  in  the  analysis  of 
test  data,  a  check  on  hypotheses,  and  hopefully  a  way  to  generate  simple  rules 
and  generalizations.  It  must  be  recognized  that  (a)  the  extent  of  the  input 


data  that  are  known  a  priori  is  limited,  since  several  key  parameters  are 
usually  not  measured,  and  (b)  the  details  of  some  of  the  important  mechanisms 
are  not  yet  well  enough  understood  to  allow  confident  predictions.  Given 
these  limitations,  one  cannot  hope  to  construct  a  truly  predictive  code. 

We  discuss  in  Section  6  optical  considerations  in  order  to  provide  a 
framework  in  which  the  results  of  the  physical  and  meteorological  modeling  can 
eventually  be  linked  to  effects  on  transmission.  Specifically,  we  discuss 
the  need  for  area,  volume,  and  settling  velocity  parameters  in  calculating 
optical  extinction  coefficients. 


2.  MODEL  ASSUMPTIONS 


Ice  particles  Introduced  in  a  supercooled  water  cloud  can  grow  rapidly  by 
vapor  diffusion  at  the  expense  of  evaporation  from  the  liquid  water.  This 
mechanism,  first  postulated  by  Findeisen,  is  fundamental  for  initiating 
precipitation,  particularly  in  the  wintertime.  It  is  due  to  the  fact  that  the 
vapor  pressure  of  liquid  water  exceeds  that  of  ice  by  as  much  as  0.2  mm  Hg  for 
a  range  of  temperatures  below  0°C,  the  difference  peaking  at  -12°C. 

The  ice  growth  rate  will  peak  at  ambient  temperatures  slightly  below  the 
peak  of  the  differential  vapor  pressure  curve  (since  the  release  of  latent 
heat  will  slightly  warm  up  the  ice  particle),  and  will  be  slower  on  either 
side  of  this  maximum.  This  growth  rate  in  turn  determines  the  shape  and  type 
of  ice  crystals  formed,  and  is  the  basis  for  ice  habit  diagrams  such  as  that 
of  Nakaya  or  Magono  and  Lee  (1966).  For  instance,  very  slow  growth  (either 
near  0°C  or  below  -25°C)  will  favor  the  thermodynamically  stable  equilibrium 
crystal  habit,  short  columns,  while  fast  growth  gives  prominence  to  kinetic 
effects  which  tend  to  accelerate  growth  of  the  prismatic  faces  and  results  in 
plate  or  branched  plate  formation  between  -10  and  -20°C.  The  degree  of 
ambient  moisture  has  a  secondary  effect  in  affecting  this  habit  selection. 
Thus,  rules  can  be  implemented  in  a  model  to  select  the  correct  crystal  habit 
depending  on  local  growth  rate. 

For  a  1  km  thick  cloud  the  temperature  difference  between  top  and  bottom 
is  between  5  and  10oC,  which  is  a  relatively  narrow  range  in  a  crystal  habit 
diagram.  Thus,  we  adopt  the  simplifying  assumption  that  only  one  crystal  type 
exists  in  the  cloud.  The  type  (plate,  needle,  or  spherical)  is  selected  a 
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priori,  based  on  Che  given  vertical  profiles.  Within  each  basic  type,  the 
degree  of  riming  by  supercooled  drops  is  calculated  in  detail,  and  it  results 
in  a  continuum  of  particle  aspect  ratios,  leading,  for  heavy  rime,  to 
spherical  graupel. 

The  origin  of  the  ice  crystals  is  one  point  of  difficulty,  since  the 
state  of  knowledge  is  still  imperfect  and  several  mechanisms  are  involved 
(seeding  from  an  overlying  "seeder  cloud",  as  originally  proposed  by  Bergeron, 
or  intra-cloud  ice  nucleation  by  contact  or  other  mechanisms).  In  any  case, 
ice  nuclei  appear  preferentially  at  the  cloud  top.  For  modeling  purposes,  we 
bypass  this  difficult  by  treating  the  nuclei  flux  at  cloud  top  in  a  parametric 
fashion. 

One  other  element  of  modeling  difficulty  concerns  the  magnitude  and 
distribution  of  the  in-cloud  updraft.  For  stratiform  or  nimbus  types  of 
clouds,  this  updraft  is  likely  to  be  in  the  range  of  0.1  to  0.5  m/s  only, 
but  since  the  fall  speed  of  ice  and  snow  particles  is  in  that  range  as  well, 
its  influence  is  essential.  Data  on  stratus  cloud  updraft  are  very  scanty, 
and  the  planned  SNOW  tests  do  not  comtemplate  its  measurement.  Therefore,  we 
treat  this  effect  parametrically  as  well,  and  will  have  to  rely  on  detailed 
comparison  of  model  results  with  those  data  which  will  be  available  (lapse 
rate,  dew  point  profile,  synoptic  weather  data,  ground  precipitation  data)  in 
order  to  construct  rules  for  selecting  the  proper  updraft  values. 

We  consider  in  the  model  a  monodisperse  water  drop  population.  The 
effect  of  polydispersity  of  liquid  water  droplets  on  their  accretion  on  ice 
particles  was  studied  by  Ryan  .  While  for  the  wide  spectrum  of  a  maritime 
cumulus  cloud,  there  was  a  factor  of  2  less  growth  if  the  mean  size  only  was 
considered,  the  effect  was  negligible  for  a  more  concentrated  continental 
cloud  spectrum. 
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Most  of  our  calculations  are  done  using  the  assumption  of  no  external  air 
entrainment  and  mixing,  although  a  mixing  length  parameter  is  specified  and  a 
sink  term  due  to  mixing  is  included  in  several  of  the  equations  (see  below). 
Aside  from  this,  the  concentration  of  liquid  water  present  is  determined  by  an 
overall  water  balance  that  includes  advection  of  vapor  and  condensation  on  ice 
particles.  The  latter  is  treated  in  a  fully  kinetic  manner,  but  the  growth  of 
individual  liquid  droplets  occurs  over  such  short  time  scale  that  its 
inclusion  as  a  kinetic  equation  would  make  the  system  of  equations  very 
stiff.  Accordingly,  above  a  threshold  of  10“  Kg/m  liquid  water  contents  it 
is  assumed  that  vapor-water  equilibrium  is  maintained.  In  the  presence  of 
large  enough  ice  concentrations,  or  if  the  supply  of  condensible  vapor  is  weak 
for  the  available  ice,  the  pressure  of  vapor  may  fall  below  water  saturation; 
this  is  detected  in  the  model  by  a  tendency  of  the  liquid  water  density,  as 
calculated  under  the  water  equilibrium  assumption,  to  drop.  When  the  10~5 
Kg/m  threshold  is  crossed,  the  assumption  of  vapor-water  saturation  is 
abandoned.  The  drops  are  assumed  to  evaporate  rapidly  (in  one  time  step)  and 
the  vapor  concentration  is  then  determined  by  the  water  balance  equation, 
without  a  liquid  term  in  it.  The  equation  that  controls  the  number  density  of 
water  drops  is,  however,  retained,  so  that  the  nuclei  at  least  are  assumed  to 
remain  present;  if  the  vapor  conservation  equation  indicates  later  a  tendency 
to  overshoot  the  water  saturation  level,  re-nucleation  is  considered  to  occur, 
and  the  excess  vapor  is  distributed  over  the  water  condensation  nuclei. 

Although  the  above  procedure  constitutes  a  satisfactory  determination  of 
the  liquid  water  contents,  the  particle  size  has  to  be  left  as  a  model 
parameter.  What  is  actually  done  is  to  specify  at  the  cloud  bottom  the 
concentration  of  cloud  condensation  nuclei.  This  serves  as  the  lower  boundary 
condition  for  a  droplet  number  conservation  equation;  the  equaton  itself 
accounts  for  time  variation,  convection  in  an  inhomogeneous  cloud,  and  loss  of 
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drops  due  to  accretion  on  ice.  The  droplet  size  then  follows  at  each  time  and 
place  from  the  water  contents  and  the  number  density. 

Accretion  (riming)  is  modeled  as  a  process  of  capture  of  droplets  by  the 
faster  falling  ice  crystals.  The  capture  efficiency  depends  quite  critically 
on  both,  drop  and  crystal  size;  we  use  curve  fits  to  the  calculated 
efficiences  of  Refs.  5  and  6- 

The  accreted  water  is  assumed  to  freeze  on  contact.  This  may  not  be  true 
at  high  accretion  rates  and  moderately  warm  temperatures,  and  some  of  the 
unfrozen  water  may  then  be  shed  off;  however,  the  details  of  this  are  not  very 
well  understood,  and  for  the  present  we  ignore  shedding.  Formation  of  ice 
splinters  upon  accretion  is  also  ignored,  as  is  crystal  breakage  by  the 
collision.  The  former  is  rarely  observed  for  droplets  smaller  than  50  pm 
diameter,  while  the  latter  requires  the  presence  of  rather  large  drops 
(diameter  >  250  pm).  Thus,  for  the  range  of  most  interest  in  winter  clouds 
(less  than  50  pm  diameter),  ice  multiplication  processes  can  be  at  least 
tentatively  neglected.  The  accreted  mass  is  assumed  to  be  added  to  the 
shorter  particle  dimension,  until  it  becomes  a  spherical  graupel  particle. 

The  density  of  accreting  particles  is  at  present  taken  to  be  0.5  g/cm8,  but  a 
more  detailed  model  is  probably  needed. 

The  primary  particle  variable  is  its  mass.  Once  this  is  known,  well 
established  dimensional  relationships  can  be  employed  to  infer  the  geometrical 
ratios  for  each  type  of  particle,  as  well  as  its  density.  The  prototype  of 
these  relationships  for  diffusional  growth  is  the  Auer-Veal  rules,  as  complied 
for  instance  by  Davis7.  Densities  are  calculated  according  to  Heymsfeld8. 

When  both  accretion  and  vapor  growth  are  present,  the  particle  diameter  is 
assumed  to  grow  purely  due  to  the  vapor  deposition,  at  a  rate  determined  by 
the  rate  of  vapor  condensation  in  the  same  manner  as  if  no  accretion  existed. 
The  accreted  mass,  on  the  other  hand  contributes  to  increase  the  thickness 
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only,  as  indicated.  This  procedure  was  also  used,  for  example,  by  Young 
(Ref .9). 


One  final  important  assumption  is  at  present  made  in  the  model,  namely, 
ice  particle  agglomeration  is  neglected.  This  is  consistent  with  our  use  of  a 
single  ice  particle  type  and  size,  since  ice-ice  agglomeration  is  driven 
mostly  by  the  difference  of  settling  velocities.  The  adequacy  of  this 
criterion  can  be  tested  by  comparing  the  ice-ice  mean  collision  frequency  to 
the  Inverse  of  the  residence  time  for  a  particle;  the  criterion  for  no 
agglomeration  is  then  found  as 


n 


i 


< 


(1) 


where  n^  is  the  ice  particle  number  density,  v^  its  mean  fall  velocity, 
with  characteristic  spread  Avj,  Aj  the  particle  frontal  area  and  H  the 
height  of  the  cloud  top.  Assuming  a  power  law  vt  *  k  ab  for  the  fall 
speed  (where  a  »  particle  radius  and  b  “  0.3-0. 6),  we  have  (A  v^/v^  ■ 
b  (A  a)/a.  For  H  »  2000m,  a  *  500  urn,  b  «  0.5  and  A  a  “  a  (as  in  an 
exponential  distribution),  the  criterion  for  no  agglomeration  gives 

nt  <  1.3  x  103  m~3  *  1.3  (liter)-1  (2) 

This  density  is  exceeded  in  glaciated  clouds  (where  however,  a  is 
probably  smaller  than  500  um)  or  in  thin  layers  of  high  concentration. 
Therefore,  for  clouds  topping  at  less  than  about  2-3  km,  agglomeration  is  not 
a  dominant  issue;  it  is  even  less  so  if  only  the  in-cloud  part  of  the  ice 
crystal  growth  is  of  interest.  As  an  example,  the  detailed  set  of  data  of 
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Ref.  10  show  no  agglomeration  in  clouds  topping  at  about  2000m,  even  when  the 
glaciated  cloud  top  reaches  ice  concentrations  of  20-40  (liter)-*  . 

However,  this  does  not  mean  that  conditions  do  not  arise  where 
agglomeration  is  important.  Cloud  top  heights  of  over  3  km  are  not  uncommon, 
and  are,  in  fact,  characteristic  of  the  most  intense  snowstorms.  In  fact,  our 
preliminary  analysis  of  the  data  of  the  SNOW  I  test  series11  indicates  the 
predominance  of  agglomeration  effects  (a  discussion  of  these  data  appear 
elsewhere  in  the  paper).  Therefore,  an  extension  of  the  model  to  cover  such 
cases  is  very  desirable  from  a  practical  viewpoint. 
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3.  FORMULATION 


The  model  consists  of  six  one-dimensional,  time-dependent  partial 
differential  -juations  for  the  following  variables:  ice  particle  number 
density  (nj),  and  mass  (m^),  accreted  mass  per  ice  particle  (ma),  ice 
particle  major  radius  (a),  droplet  number  density  (n,j)  and  droplet  mass 
(m^)*  If  mixing  with  outer  air  is  neglected  (case  of  very  extensive  stratus 
layers),  the  equations  are: 
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3a  .  3a 
t: —  +  v.  r — 
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Here  vj  and  vd  are  the  net  (upwards)  velocity  of  ice  particles  and 
droplets  respectively: 
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v  —  V 

a  i,  seed 

V  —  v 

a  d,  settl 


(9) 

(10) 


where  Va  is  the  prescribed  updraft  velocity,  and  vj>setti,  vd>setti  are 
the  ice  and  droplet  settling  velocities,  respectively.  The  vapor  density  pv 
is  taken  as  the  water-equilibrium  density  (from  a  curve  fit  to  the 
Clausius-Clapeyron  equation)  whenever  ndmd  >  10  Kg/m  .  Since  pv  then 
depends  only  on  the  ambient  temperature  (whose  profile  is  prescribed),  Eq .  (6) 
can  be  viewed  as  governing  the  droplet  mass  ny .  On  the  other  hand,  if 
ndmd  <  10  Kg/m  ,  ny  is  set  equal  to  zero,  and  Eq .  (6)  is  turned  around 
to  become  a  conservation  law  for  pv  (breaking  the  connection  of  pv  to  T) . 
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The  diffusive  rate  of  growth  of  ice  particles  is  given  by 


dm . 

^dt^diff 
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where  C  is  the  electrostatic  capacitance  of  the  particle  (a  function  of  its 
size  and  shape),  e0  *  8.85  x  10~12  in  MKS  units  is  the  permittivity  of 
free  space,  SV£  is  the  supersaturation  relative  to  ice,  R-8.314  Joule/(Mole 
°K)  is  the  gas  constant,  esatji  is  the  ice  saturation  vapor  pressure,  M,, 
is  the  molecular  weight  of  water,  Ls  the  heat  of  ice  sublimation  and  Dv' , 

Ka'  are  effective  vapor  diffusivity  and  air  thermal  conductivity, 
respectively.  These  effective  transport  properties  incorporate  correction 
factors  for  both  free-molecular  flow  effects  and  high  Reynolds  number  effects; 
details  can  be  found  for  example  in  Ref.  (12),  Ch.  13.  Equation  (11)  accounts 
in  its  formulation  for  the  reduction  of  condensation  rate  due  to  particle 
heating  by  the  release  of  latent  heat.  The  ice  supersaturation  SVfi  is 
calculated  by 
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where 

e  =  p  R  T  (13) 

00  V  V 

and  e  stands  for  the  prevailing  pressure  of  vapor  (e»  <  esat>w  (Ta))« 
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Accretion  growth  for  plate-like  particles  is  given  by: 
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where  Ec  is  the  collection  efficiency,^  a<j  and  aj  are  the  drop  and 
plate  radii,  respecf  iv-  :.y,  and  is  the  local  water  content  of  the  cloud. 
Similar  expressions  are  used  for  other  particle  shapes. 

Particle  1  droplet  fall  velocities  are  calculated  using  standard 
literature  value*  of  the  drag  coefficient  for  the  various  shapes  and  sizes. 
The  calculation  t^n.iea  the  form 
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can  be  calculated  directly  in  terms  of  the  particle  mass  inj,  f  is  a  known 
function,  depending  on  geometry  only,  and 


o  v  d 
w  a  FALL 


Re 


(17) 


is  the  particle  Reynolds  number,  from  which  the  fall  velocity  follows. 
Functions  f  for  both  oblate  and  prolate  spheroids  (as  well  as  for  spheres) 
have  been  incorporated  in  the  model. 
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Equation  (8)  expresses  the  assumption  that  only  vapor  condensation  can 
increase  the  particle  major  dimension,  a.  The  quantity  (dm^/da^tf f  is 


calculated  from  the 

expressions  in  Refs. 
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and 
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Mixing  with  external  air  is  accounted  for  through  a  simple  device.  It  is 
assumed  that  mixing  is  proportional  to  updraft  strength  va,  inversely 
proportional  to  cloud  horizontal  extent,  and  proportional  to  the  difference  in 
densities  of  the  mixing  quantities  (ice,  water,  and  vapor  densities).  No 
enthalopy  mixing  is  explicitly  included,  since  the  temperature  profile  is 
directly  an  input  to  the  model.  Introducing  a  “mixing  length”  £m  (of  the 
order  of  5  times  the  horizontal  radius  of  the  cloud,  according  to  Ref.  13), 
the  additional  terms  on  the  right-hand  sides  of  Eqs.  (3), (4),  and  (6)  are 
respectively 


P  ) 
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(3a) 
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(6a) 


It  can  be  seen  by  comparison  of  these  terms  with  the  convective 
derivative  terms  in  the  corresponding  equations  that  mixing  effects  become 
significant  only  for  clouds  with  thickness  comparable  to  their  diameter;  for 
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extensive  stratus  (£m  »  cloud  thickness),  mixing  is  only  significant  near 
the  cloud  edges. 

The  "natural"  boundary  conditions  of  the  problem  are:  ice  flux  or  ice 
number  density  and  ice  particle  mass  at  cloud  top,  water  droplet  number 
density  and  zero  droplet  diameter  at  cloud  base,  zero  accreted  mass  and 
prescribed  particle  diameter  at  cloud  top,  and  water-saturation  vapor  density 
at  cloud  base. 

In  addition,  the  numerical  procedure  requires  the  remaining  boundary 
conditions  to  be  specified  as  well  (on  ice  mass,  diameter  and  accretted  mass 
at  cloud  base,  droplet  density  and  diameter  and  vapor  density  at  cloud  top). 
We  choose  zero  flux  gradient  conditions  for  these  quantities,  indicative  of 
the  fact  that  no  changes  (or  only  very  slow  ones)  are  expected  outside  the 
cloud . 

The  initial  conditions  are  normally  chosen  to  correspond  to  an  ice-free 
supercooled  cloud  in  a  nonprecipitating  steady  state.  Some  inconsistency  is 
incurred  in  that  precipitation  would  be  required  to  balance  the  vapor 
advection  by  the  updraft,  which  is  assumed  constant  in  time.  This  initial 
inconsistency  disappears  gradually  as  ice  phase  precipitation  develops  from 
the  cloud  top  down. 

The  system  of  partial  differential  equations  is  solved  in  time  and  space 
by  the  McCormick  time-split  method,  accurate  to  first  order  in  time  and  to 
second  order  in  space.  Typical  mesh  sizes  are  20m  and  5  sec  (if  droplet 
growth  were  modeled  kinetically,  a  much  shorter  time  step  would  be  required) . 
A  slow  growing  numerical  instability  of  short  wavelength  (two  grid  spaces)  is 
noticed  for  some  conditions.  A  periodic  smoothing  of  the  profiles  using  a 
local  filter  has  been  found  to  correct  this  problem. 
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4.  APPLICATIONS  OF  THE  MODEL 


We  present  in  this  section  some  examples  of  the  sort  of  results 
obtainable  with  the  model  just  described.  Although  comparison  to  data  is  made 
difficult  by  the  lack  of  detailed  space-time  measurements  of  in-cloud 
developments,  some  qualitative  comparisons  will  be  made  to  the  data  of  Ref. 

10.  as  well  as  (in  the  following  section)  to  SNOW  I  test  data. 

Figures  4  through  7  present  in  detail  the  calculated  time  history  of  a 
snowstorm  where  most  of  the  effects  included  in  the  model  can  be  observed. 

The  input  data  are  as  follows: 

Cloud  thickness:  From  1000  to  2000m  in  altitude 

Temperature  profile:  Linearly  decreasing  from  -11°C  to  -17°C 
(8°C/km).  This  is  close  to  moist  adiabatic  conditions. 

Updraft:  17cm/sec,  tapering  to  8.5  cm/ sec  in  the  last  60m. 

Ice  nuclei  flux:  1000/(m  sec)  at  cloud  top 

Ice  nuclei  size:  200  y  radius  at  cloud  top 

8  3 

Droplet  density:  2  x  10  m“  at  cloud  base 

Mixing  length:  3000m 

Ice  particle  type:  plates  (consistent  with  growth  at  about  -15°C) 

Some  of  these  data  were  selected  to  match  the  conditions  reported  in 
Ref.  10  (temperatures,  altitudes,  droplet  size,  crystal  type).  In  that  work, 
clouds  in  different  stages  of  development  were  observed  during  a  two  day 
period  and  in-cloud  data  were  sampled  (ice  and  water  contents,  crystal 
concentration  and  type).  Although  the  clouds  are  classified  as  "incipient”. 
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"mature"  or  "decaying",  no  information  can  be  obtained  regarding  their  earlier 
history  or  other  conditions  (updraft,  nuclei  flux). 

Figure  1  shows  the  development  of  the  ice  density  profile.  After  1000 
s  from  ice  seeding  initiation,  the  ice  has  penetrated  some  800ra  into  the 
cloud;  accounting  for  the  updrafts,  this  represents  a  settling  velocity  of 
some  95  cm/sec,  which  is  typical  of  rimed  sub-millimeter  plates.  The 
development  from  this  time  on  is  quite  regular,  with  a  gradual  increase  in  ice 
concentration,  approaching  about  3  particles/lirer .  This  increase  indicates 
accumulation  due  to  slower  fall  velocities  as  time  increases;  we  will  later 
trace  it  to  the  gradual  disappearance  of  rime  on  the  particles,  as  the  cloud 
glaciates. 

Figure  2  shows  the  development  of  ice  crystal  radius.  The  gradual 
increase  with  time  corresponds  also  to  the  slower  fall  velocity,  which  gives 
each  particle  more  residence  time,  and  allows  more  growth  by  vapor  diffusion. 
By  about  11,000  sec,  an  acceleration  of  this  growth  is  discernible;  this  can 
be  correlated  with  the  accumulation  of  crystals  at  cloud  top,  visible  in  Fig. 

1  at  11,000  sec,  and,  as  will  be  seen  later,  with  the  arrival  of  glaciation 
there.  The  disappearance  of  accretion  leads  to  this  local  accumulation  due  to 
the  very  slow  fall  of  the  small  ice  crystals  present  at  the  top. 

Figure  3  shows  the  accreted  mass  per  particle.  The  main  accretion  event 
occurs  between  1000  and  3000  sec,  and  the  1000  s  profile  shows  a  significant 
amount  of  riming;  at  the  peak  of  this  profile,  about  two  thirds  of  a  particle 
mass  is  in  the  accreted  droplets.  As  time  progresses,  the  accreted  mass 
gradually  disappears,  and  is  insignificant  at  11,000  sec. 
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ALTITUDE  (km) 


ICE  NUMBER  DENSITY  (PER  LITER) 


Fig-  1 


Ice  Number  Density,  Ice  Particle  Flux  at  Cloud  Top  is  1000  m~2s_1 , 
Water  Droplet  Density  at  Cloud  Bottom  is  200  cm-^ . 
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ALTITUDE  (km) 
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Fig.  3  Af2crf‘/ed  Mass  Per  Particle.  Ice  Particle  Flux  at  Top  is  1000 
s~  .  Water  Droplet  Density  at  Cloud  Bottom  is  200  cm-3 


27 


The  corresponding  developments  regarding  the  liquid  phase  are  presented 
in  Figs.  4,  5  and  6.  In  Fig.  4  we  can  observe  the  strong  depletion  of  water 
droplets  by  the  process  of  riming  at  t  *  3000-5000  sec,  particularly  in  the 
upper  third  of  the  cloud,  where  the  drops  have  diameters  of  about  20  Um,  and 
are  at  the  peak  of  their  accretion  efficiency.  The  reduction  of  their 
efficiency  down  the  cloud,  as  the  droplets  become  smaller  (see  Fig.  5)  can  be 
easily  traced;  it  overcomes  a  partially  compensatory  effect  associated  with 
the  higher  capture  efficiency  of  the  growing  ice  crystals.  As  glaciation 
proceeds  by  evaporation  of  the  droplets,  the  number  density  profile  recovers 
to  its  original  level,  although  by  now  it  represents  merely  "dried  out" 
nuclei.  Figure  5  shows  the  slow  progress  of  glaciation.  The  droplets  first 
evaporate  at  cloud  base,  where  the  most  ice  surface  is  available;  however, 
this  is  sensitive  to  parameters  such  as  lapse  rate  and,  possibly, 
agglomeration  (not  treated  in  the  model).  The  survival  of  liquid  water  at 
cloud  top  is  an  artifact  of  the  model,  due  to  fictitious  accumulation  as  the 
updraft  tapers  off.  Of  course,  lateral  divergence  must  exist,  to  compensate 
for  this  slowdown,  and  this  should  introduce  an  equivalent  one-dimensional 
"sink"  in  the  conservation  equations.  The  total  liquid  water  contents  of  the 
cloud  is  shown  in  Fig.  6,  where  a  progression  from  a  profile  that  reaches 
”  1  g/m  at  cloud  top  to  essentially  zero,  can  be  seen  to  occur. 

The  vapor  density  is,  of  course,  a  driving  factor  throughout  these 
developments.  Its  evolution  is  seen  in  Fig.  7.  The  initial  profile  is  simply 
a  reflection  of  the  temperature  profile  assumed,  since  water  saturation  exists 
initially.  At  the  latter  times,  when  the  cloud  is  glaciated,  the  vapor 
density  falls  very  slightly  below  water  saturation.  A  lower  bound  would 
probably  exist  at  the  ice  saturation  level  (85-90%  of  water  saturation  at 
these  temperatures),  but  this  limit  is  not  even  approached;  this  means  enough 
supersaturation  still  exists  relative  to  the  ice  phase  to  ensure  continued 
growth  (see  Fig.  2)  in  these  later  stages.  In  turn,  this  prevents 
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DROPLET  RADIUS  (km) 


Droplet  Radius.  Ice  Particle  Flux  at  Cloud  Top  is  1000  m' 
Water  Droplet  Density  at  Cloud  Bottom  is  200  cm-3 


Fig.  6 


-  Liquid  Water  Content.  Ice  Particle  Flux  at  Top  is  1000  m-2s_1. 
Water  Droplet  Density  at  Cloud  Bottom  is  200  cm-3 
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substantially  greater  ice  accumulation  at  the  cloud  top,  which  would  occur  if 
slower  growth  happened  there.  This  kind  of  accumulation  is  seen  in  the  data 
of  Ref.  10  for  the  mature  and  decaying  stages;  the  authors  of  Ref.  10 
attribute  it  to  some  sort  of  ice  multiplication  process.  We  believe  it 
possible  that  large-scale  accumulation  may  also  result  from  slow  growth,  if 
conditions  leading  to  cummulative  vapor  depletion  are  established.  Continued 
work  in  improved  modeling  of  the  cloud  top  region  is  needed  to  clarify  this 
point . 

A  pre-condition  for  glaciation  is  a  sufficiently  high  ice  particle 
concentration,  or,  conversely,  a  sufficiently  low  lapse  rate,  so  that  the  ice 
phase  may  absorb  all  of  the  excess  vapor  carried  by  the  updraft.  Criteria  for 
this  threshold  have  been  established,  for  example  by  Jiusto  (Ref.  14);  for  our 
conditions  it  seems  to  correspond  to  about  3  particles/liter,  as  in  the 
previous  example. 

To  verify  this  limit,  we  repeated  the  calculations  for  similar 
conditions,  except  that  the  ice  flux  at  the  top  was  reduced  by  a  factor  of  2 
(to  500  m  see*  ).  Figure  8  shows  the  development  of  the  liquid  water 
contents  with  time  in  their  case;  there  is  a  strong  reduction,  traceable 
mostly  to  the  continued  accretion  of  droplets  on  ice  particles,  but  the 
droplets  show  no  tendency  to  evaporate  in  the  presence  of  about  1  ice 
crystal/liter.  Notice  the  ice  concentration  has  decreased  by  a  factor  of  3, 
while  the  flux  was  reduced  only  by  2;  this  is  due  to  the  faster  fall  speed 
(0.8-1  m/sec  vs  0.5-0. 8  m/sec)  of  the  particles  in  the  second  case,  which  in 
turn  results  from  the  continuation  of  riming  throughout  the  cloud 
development.  Figure  9  illustrates  this;  about  half  the  mass  per  particle  is 
rime,  even  after  11,000  sec.  Also  associated  with  this  is  the  smaller  radius 
of  the  ice  particles  (0.8  pm  at  cloud  base),  due  to  their  smaller  residence 
time.  Altogether,  then,  the  riming  behavior  of  the  ice  particles  cooperates 
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LIQUID  WATER  CONTENTS  (g/V) 


Liquid  Water  Content.  Ice  Particle  Flux  at  Top  is  500  m_2s' 
ater  Droplet  Density  at  Cloud  Bottom  is  200  cm-3 


ALTITUDE 


Fig.  9  -  Accreted  Mass  per  Particle.  Ice  Particle  Flux  at  Top  is  500  m” 
Water  Droplet  Density  at  Cloud  Bottom  is  1000  cm- 
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in  making  a  fairly  fast  transition  from  glaciating  to  non-glaciating 
conditions  as  the  ice  flux  is  reduced. 

Going  to  the  opposite  end  of  the  concentration  scale,  we  tested  a  case 

p 

with  2000  crystals/(m  sec)  (twice  the  flux  of  the  basic  case). 

Simultaneously,  the  diameter  of  the  droplets  was  reduced  to  about  12  um,  to 
reduce  accretion  effects.  This  resulted  in  a  very  rapid  glaciation  of  the 
whole  cloud  between  t  =  3000  sec  (when  the,  now  lightly  accreting  crystals, 
first  reach  the  cloud  base)  and  about  t  =  4000  sec.  The  vapor  density 
profiles  for  this  case  are  shown  in  Fig.  10.  The  vapor  is  clearly  below  water 
saturation,  but  still  only  about  half  way  to  ice  saturation.  It  would  appear 
that  in  order  to  prevent  vapor  diffusion  growth  in  the  upper  part  of  the 
cloud,  a  crystal  seeding  rate  even  two  to  three  times  higher  is  required; 
under  such  conditions,  accumulation  at  the  top  can  be  expected. 
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5.  OBSERVATIONS  ON  THE  SNOW  I  DATA 

Due  to  the  apparent  existence  of  agglomeration  in  at  least  parts  of  the 
snowstorms  of  January  and  February  1981,  as  reported  in  Ref.  (11),  we  have  not 
attempted  a  direct  simulation  using  the  numerical  code.  We  will  limit 
ourselves  for  now  to  a  few  comments  on  the  implications  of  the  data. 

Figure  11  shows  vertical  profiles  taken  shortly  after  the  start  of  snow 
on  January  22.  There  is  a  nearly  isothermal,  relatively  thin  cloud  which 
appears  to  be  ice,  rather  than  water  saturated.  Glaciation  is,  indeed, 
favored  by  low  lapse  rates,  because  a  moderate  ice  particle  concentration  can 
then  absorb  the  small  amount  of  vapor  made  available  by  the  updraft.  However, 
the  clouds  may  be  reverting  to  liquid,  since  light  to  medium  riming  is 
reported  later  on  the  day.  At  the  early  times,  the  cloud  could  be  classed  as 
a  "seeder-only"  system,  in  the  nomenclature  of  Jiustr-  and  Weickmann 
(Ref.  (15)). 

Figure  12  corresponds  to  the  peak  of  the  snowstorm  of  February  8-9, 

1981.  The  in-cloud  temperature  follows  closely  a  moist  adiabatic  line, 
showing  a  well  developed  stratus  system.  The  cloud  is  quite  thick,  water 
saturation  being  maintained  to  at  least  3200m,  with  a  low  cloud  base 
(*400m).  The  temperature  range  goes  from  -5°C  to  -25°C,  which  explains  the 
wide  variety  of  snow  particles  observed  at  the  time.  The  very  small  columns 
reported  may  actually  have  formed  near  cloud  base,  in  the  -5  to  -10°C  range, 
while  the  more  complex  and  large  crystals  must  come  from  (or  at  least  have 
gone  through)  the  mid-regions  of  the  cloud,  around  the  maximum  growth  rate 
range  (-12  to  -18°C) . 
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Fig.  11  -  Data  From  Airsonde  Flight,  1130,  22  January  1981  (From  Reference 
(11),  With  Moist  Adiabats  Added) 
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Fig-  12  -  Data  From  Airsonde  Flight,  2210,  8  February  1981  (From  Reference 
(11),  With  Moist  Adiabats  Added) 
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The  data  of  snow  concentration  and  precipitation  rate  are  shown  in 
Fig.  13  plotted  versus  each  other  for  comparison  with  the  empirical 
correlation  of  Sekhon  and  Srivastava  (Ref.  (16)).  Most  points  cluster  about 
the  correlation  line,  with  corresponding  fall  velocities  from  1  to  1.4  m/'sec. 
For  those  points  deviating  clearly  to  the  right  (high  fall  velocities,  2-5 
m/sec),  one  can  recognize  conditions  of  heavy  riming  in  the  descriptions 
reported . 

It  is  of  interest  that  Passarelli  (Ref.  (17))  has  been  able  to  obtain  a 
theoretical  prediction  of  the  Sekhon-Srivastava  line  based  on  the  existence  of 
a  limiting  “equilibrium"  growth  regime,  in  which  both  agglomeration  and  vapor 
growth  (with  their  opposite  effects  on  polydispersity)  are  acting 
simultaneously.  Based  on  the  equations  of  Ref.  (17)  the  equilibrium  growth 
regime  is  only  reached  after  a  fall  distance  of  the  order  of  £Vp/Kl ,  where 
i,Vp  is  the  length  seal’  for  variation  of  vapor  pressure  in  the  cloud  and  Kj 
is  only  a  function  of  the  exponents  for  the  fall  speed-diameter  and  vapor 
growth-diameter  relationships;  typical  values  of  K;  are  about  0.3.  Using  the 
Claus ius-Clapeyron  equation,  then,  we  obtain  a  "development  distance" 


(19) 


where  Rv  is  the  vapor  gas  contant,  Ls  the  heat  of  sublimation  and  T  the 
lapse  rate.  For  T  «  -16°C,  T  »  7°C/Km,  one  finds  i  -  5000m.  Thus,  only  for 
rather  thick  strata  can  one  expect  strict  equilibrium  growth,  although  a  good 
approximation  seems  to  be  obtained  in  only  a  fraction  of  the  above  distance. 
In  any  event,  Fig.  13  seems  to  confirm  the  importance  of  agglomeration  for  the 
case  of  the  February  8-9  storm. 
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Fig .  13  -  Triangles  are  Feb.  8-9,  1981  Storm  Data  (Time  of  Day  Indicated) 
Solid  Dot  is  From  Jan.  22,  1981  Storm  The  Solid  Line  is  Sekhon- 
Srivastava  Correlation  (Ref.  16)  Dotted  Lines  are  Model  Results 

(A)  Case  500/2/8  (B)  Case  500/5/8  (C)  Case  1000/2/8 

(D)  Case  1000/1/8  (E)  Case  2000/10/8 
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Also  shown  in  Fig.  13  are  the  results  of  our  model  calculations  at  cloud 

base,  for  five  sets  of  conditions.  These  are  denoted  by  the  code  a/b/c,  where 

2  8 

a  is  the  crystal  flux  (per  m  ,  per  sec),  b  is  10”  x  the  droplet  density  at 
cloud  base  (m~  ),  and  c  is  the  lapse  rate  in  the  upper  400m  of  the  cloud,  in 
°C/km  (the  lower  600m  are  in  all  cases  taken  at  8°C/km).  The  trends  predicted 
are  correct,  but  our  calculated  fall  velocities  tend  to  be  lower  than  observed 
in  the  SNOW  data  and  than  implied  by  the  correlation  data.  This  is  most 
likely  a  result  of  the  neglect  of  agglomeration,  which  increases  particle  mass 
somewhat  faster  than  exposed  area. 
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6.  SNOW  OPTICS  TRANSMITTANCE  CALCULATIONS 


In  this  section  we  develop  a  relationship  between  extinction  coefficients 
in  the  visible/infrared  to  snowflake  physical  characteristics.  We  begin  with 
the  equation  for  the  transmission  through  falling  snow 


t  =  exp  (-6L) 


(20) 


where 

L  =  path  length 

d  =  extinction  coefficient  *  / A  N(A)Qext  (nice»A,  A)  dA 
N(A)  =  number  of  snowflakes  with  area  A 
Qext  =  extinction  efficiency 
nice  =  refractive  index  of  ice 
A  =  wavelength 

The  equation  for  conservation  of  mass  gives 

M  =  M 

water  ice 


where 


UU 


Mvacer  ”  melted-ice  water  equivalent  of  snow 


Mfce  “  mass  of  ice 


P  R 

water 


Pice  /VN(VMV>  dv 


(21) 


pwater  ■  density  of  water 

R  =  precipitation  rate 

p ice  ”  density  of  ice 

N(V)  »  number  of  snowflakes  with  volume  V 

v(V)  *  settling  velocity  of  snowflake  with  volume  V  . 


Dividing  Eq .  (21)  by  Eq .  (20)  and  solving  for  8: 
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water 
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(22) 


For  a  monodisperse  particle  distribution,  N(V)  *  N06(V-V0)  and  N(A)  « 
Nod(A-A0),  in  which  case  Eq.  (22)  becomes: 
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(23) 


For  visible/infrared,  Qext  3  2,  so  that  a  relatively  simple  formula  for  B 
results : 


•w 


tJ 


water 


ice 


(24) 


Data  exists  (e.g.,  Refs.  2,18,19)  on  values  of  8,  A,  V,  and  v  for 
different  snow  types.  The  data  points  shown  in  Figures  14  and  15  were 
digitized  for  subsequent  curve  fitting.  (The  original  curves  shown  in  Figure 
14  and  15  are  drawn  with  unknown  curve  fits.) 

Appendix  A  shows  the  curve  fit  results  for  fall  velocity  v  versus 
dimension  for  different  snow  crystal  type.  Appendix  B  gives  the  results  for 
mass  M  (actually  the  3/2  power  of  the  linear  dimension  designated  "r"  in 
Figure  15). 
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Fig.  14 


d  *  Dimension  of  crystal 


Snowflake  Fall  Velocity  Data  Versus  Dimension  for  Different  Snow 
Crystals.  Original  Data  with  Unknown  Curve  Fits  (Ref. 2). 


A 


Snowflake  Mass  Versus  Dimension.  Original  Data  With  Unknown  Curve 
Fits  (Ref.  2) 


7.  CONCLUSIONS  AND  RECOMMENDATIONS 


A  relatively  simple  model  has  been  developed  that  can  be  used  for 
exploration  of  the  intricate  sensitivities  and  interrelationships  affecting 
the  nature  and  intensity  of  snowfall.  Although  the  model  can  be  of  use  in  its 
present  form,  preliminary  evaluations  against  data  indicate  the  need  for 
several  improvements  and  extensions. 

(a)  Inclusion  of  ice-ice  agglomeration 

(b)  Elimination  of  fictitious  accumulations  due  to  updraft  dropoff. 

(c)  Inclusion  of  wet  snow  effects 

The  work  has  also  led  to  the  awareness  that  more  detailed  in-cloud 
measurements  would  greatly  enhance  our  understanding  and  predictive 
capabilities.  Models  like  the  present  one  could  be  used  to  evaluate  and 
correlate  these  detailed  data. 

Finally,  a  simple  relationship  between  extinction  in  the  visible/infrared 
and  snowflake  physical  characteristics  provides  a  framework  in  which  the 
snowflakes  resulting  from  calculations  performed  by  the  model  described  above 
can  be  related  to  optical  transmission  data. 
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APPENDIX  A 


DIMENSION  OF  CRYSTAL 


VS. 


FALLING  VELOCITY 
FOR  SIX  TYPES  OF  SNOW 


Figure  A-l .  Crystal  dimension  (mm)  vs  falling  velocity  (cm/s)  for  needles 


FALLING  VELOCITY 


DIM  VS  FALLING  VELOCITY  #2  PLANE  DENORITIC 


DIMENSION  OF  CRYSTAL 


Figure  A-2 .  Crystal  dimension  (mm)  vs  falling  velocity  (cm/s)  for  plane 
dendritic. 
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FRLLING  VELOCITY 


DIM  VS  FRLLING  VELOCITY  #3  SPATIAL  DENDRITIC 


DIMENSION  OF  CRYSTAL 


Figure  A-3.  Crystal  dimension  (mm)  vs  falling  velocity  (cm/s)  for  spatial 
dendritic. 


FALLING  VELOCITY 


DIM  VS  FALLING  VELOCITY 
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Figure  A-5.  Crystal  dimension  (mm)  vs  falling  velocity  (cm/s)  for  crystal 
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APPENDIX  B 


DIMENSION  OF  CRYSTAL 
VS. 

R**  (3/2) 

FOR  SIX  TYPES  OF  SNOW 
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*»  (3/2) 


DIM  VS  R**  (3/2)  #1  NEEDLE 


cr 


Figure  B 


-1.  Crystal  dimension  (mm)  vs  r  (mm)  to  the  3/2  power  for  needles. 
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DIMENSION  OF  CRYSTAL 


Figure  B-2.  Crystal  dimension  (mm)  vs  r  (mm)  to  the  3/2  power  for  plane 
dendritic. 
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OIM  VS  R*»  £3/2)  #3  SPATIAL  OENDRITIi 


DIMENSION  OF  CRYSTAL 


Figure  B-3.  Crystal  dimension  (mm)  vs  r  (mm)  to  the  3/2  power  Eor  spatial 
dendritic. 
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R*«  C3/2) 


DIM  VS  R*« 13/23  #5  CRYSTAL  WITH  DROPLETS 


DIMENSION  OF  CRYSTAL. 


Figure  B-5.  Crystal  dimension  (mm)  vs  r  (mm)  to  the  3/2  power  for  crystal 
with  droplet. 
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